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ABSTRACT 


In this article the formal process for determining the higher order 
perturbations of the orbital elements is developed by using the methods 
of Krylov- Bogoliubov and Poincare". Such a development is necessary, 
for example, in the lunar problem where very high order perturbations 
have to be determined. The differential equations are formed for the 
elements which are affected only by the long period and the secular 
terms. The problem of determining these elements, as well as elimi- 
nating the short period effects, is reduced to solving a set of partial 
differential equations, step by step. By developing the displacement 
operator into a series of the differential operators of Faa de Bruno we 
can write these equations in a concise form. 
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ON THE HIGH ORDER EFFECTS IN THE METHODS OF 
KRYLOV-BOGOLIUBOV AND POINCARE 

by 

Peter Musen 

Goddard Space Flight Center 


INTRODUCTION 

In this article we develop the formalism for determining the general perturbations of higher 
orders in celestial mechanics by the methods of Krylov- Bogoliubov (1961) and of Poincare (1892) 
and von Zeipel (1916). The solution is obtained in terms of the Krylov- Bogoliubov averaging op- 
erator, of Faa de Bruno (1855) differential operators, and of the integrating operator. 

In the method of Krylov- Bogoliubov, as in that of Poincare, the final goal is to eliminate the 
short period effects and derive the elements affected only by the long period and secular perturba- 
tions. The original work of Krylov and Bogoliubov was influenced by the problems of celestial 
mechanics. A close look at Le Verrier's (1856) method for the general perturbations, reveals the 
same basic idea, but the method of Krylov and Bogoliubov achieved fame under its present name 
because of its extensive application to other problems of theoretical physics. In most such prob- 
lems there is no need to compute the effects of higher orders: Normally only the effects of the 
first and second orders, rarely of the third order, are computed. The standard presentation of the 
method does not go beyond these limits. However, this accuracy is insufficient from the standpoint 
of celestial mechanics. In the lunar problem, to secure the necessary accuracy of the long period 
terms, we must go up to the ninth order with respect to the ratio of mean motions of the satellite 
and of the sun. 

Thus, the formalism of the Krylov- Bogoliubov method must be extended to cover such cases 
and especially to provide for the determination of the long period effects of higher orders. The 
long period and the secular effects are chiefly responsible for the behavior and stability of the 
orbits of the celestial body, and their accurate determination is of great importance. The positive 
characteristic of the method of Krylov and Bogoliubov is that the canonical form of the equations of 
motion is not required, and thus the method can be applied to a much wider range of problems than 
the method of Poincare. However, the number of partial differential equations to be solved in the 
process of eliminating the short period terms increases, as compared to the method of Poincare. 

It is the price paid for extending the domain of applicability. 

In the method of Poincare the equations of motion have the canonical form, and the problem of 
eliminating the short period terms from the coordinates and momenta reduces to eliminating such 
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terms from the Hamiltonian by means of an appropriate canonical transformation. Assuming that 
the characteristic function S of this transformation is developable into a power series with re- 
spect to a small parameter, we reduce the determination of S to the solution of a chain of partial 
differential equations, step by step. Recently, Giacaglia (1964) has established the general form 
of these equations. 

We show in this work that the partial differential equations of the method of Poincare 7 take an 
especially concise form if written in terms of Faa de Bruno (1855) operators. 


HIGHER ORDER PER PERTURBATIONS IN THE KRYLOV-BOGOLIUBOV METHOD 

Consider the system of vectorial differential equations 


dx 

dt 


X(x; y, r\) , 


( 1 ) 


dy 

d t 


A (x ) + Y(x ; y , >i), 


(2) 


dt 


H(x; y , n) , 


(3) 


where x, Y, H are periodic in vectors y and »i with the period 2 n in each component. These 
vectors are assumed to be developable in powers of a small parameter. We have 


x = X j y ’ ’ll. 

j=i 


Y = 2^ Y. (x; y, n ), 


H = ^ ' (C (x; y, t| ") , 

i=l 

where the functions X. , Y , H are of the form 

j ’ i ’ i 

F = F n _ v (x) exp i (n -y + v -t]), 

n ,v 

where n and v are vectors whose components are integers. 


(4) 


(5) 


(6) 


(7) 
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The terms in (7) are: 


the short periodic if n ^ 0, 

the long periodic if n=0 but « / 0; 

the secular if n = 0 and v = 0. 

The averaging operator M extracts the long period and the secular terms from (7). Thus 


MF 


L 


F 0v (x) exp i vt| 


In addition to the Krylov-Bogoliubov operator M, it is also convenient to use the operator P which 
extracts the short period terms only: 


PF = / F n v (x) exp i (n-y + vt]) 

n^O , v 

In the further exposition we make use of the partial del-operators 3/3x , B/3y, 3/Bn and introduce 
a partial differential equation of the form 



PF. 


Evidently 





n/0,v 


exp i (n • y + v -q) 


Introducing the integrating operator Q we can write 


v = QPF. 


Let us determine the transformation 


+ a(x*; 

y*. 

(8) 

+ b (x *; 

y*. n*). 

(9) 

+ P(x*; 

y*. n*) 

(10) 


in such a way that the differential equations for the new variables 
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dx * 

dt 


X* (x*; n*)> 


(ID 


dy * 

dt 


X (x *) + Y*(x*; n*). 


( 12 ) 


drf 

dt 


H* (x*; t,*) 


(13) 


do not contain the new short period argument y ♦ ; we put a dash in place of y * to emphasize its 
absence. 


We shall determine the formal developments 


a 



H p = L ».• 

;=i 


X* 



X*, 




(14) 


(15) 


in such a way that (11) - (13) have the prescribed form. It follows from these equations that the 
operator d/dt can be written in the form 


where 


D = X*-- 


3x* 


+ Y* 


dy * 


H* 


From (8) - (10) and (11) - (13) we have 


= X* + A-^- + D». 
dt 3y* 


(16) 


= * + Y* + A- — + Db, 

d t 3y * 


(17) 


dq 

dt 


H* + A-^L + DjJ. 

By * 


(18) 
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Introducing the displacement operator 


1 + T(x*; y*. i]*) = exp 


/ d , d „ d \ 

a- + b- + p- 

\ dx* dy* Bt] */ 


We can write (1) - (3) as 


^ = (1 + T) X(x*; y *. n*), 

d t 


( 19 ) 


dy _ 


Y - = (1 + T) [A (x * ) + Y(x*; y*. r,*)] , 


( 20 ) 


P- = (1 + T) H(x*; y*, n *). 
d t 


( 21 ) 


Comparing (16) - (18) with (19) - (21) and changing the notation, we then have 


A • = (X - X*) + (TX - Do). 

dy 


A-^- = (Y-Y*) + TA + (TY-Db), 

dy 


A-i = (H - H*) + (TH - Dp). 

dy 


( 22 ) 

(23) 

(24) 


Making use of (14), we can represent 1 + T in the form 


1 + T = exp 


L - £ v 


J=0 


(25) 


where we put 


_d_ 

’ dx 


_d_ 

i dy 


^ i u _ d 

0i ‘^ + b ‘'“ + 


The operators T. are polynomials in 5j, S 2 , ••*. They can be decomposed into the sums 


r , = £ T <-> 


where T are homogeneous and of degree k with respect to the s-operators. Making use of the 
expressions obtained by Faa de Bruno (1855) for the higher derivatives of a function depending upon 
another function, we obtain 
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i,i 


2.1 


3. 1 


3,2 


4.1 


4 , 2 


5 , 1 


5,3 


5,5 


6,1 


6,2 


6,3 


6,4 


Is? 

2 1 


m 2 


i 8* 
6 1 


+ 


I 8 ? 8 * 


— sf 

24 1 

S, 


S l 5 4 + 


1 8 2 S 

2 1 3 


l *?*. 


— — - sf 

120 1 


S ! S 5 + 


- S ? S 4 
2 1 4 


1 84 83 


1 

2 


8 


2 

2 


8 2 S 3 


+ 7*1*5 


S 2 S 4 4 5 83 


+ 828283 4 4 


4- s? s* 


S 


3 

2 


6 



6,5 


6,6 


— S? 8, 

24 1 2 


— 8® 
720 1 


T 7,l = S 7 


7,2 


7,3 


7,4 


* 1*6 + S 2 S 5 + S 3 S 4 


7,7 


I 8 2 S 5 + 8j S 2 8 4 + | S t 8| + I S 2 S 3 


I S3 s 4 + i 8 2 s 2 S 3 + I 8, S3 


£ 8 ? 8 > ^ £*2*5 


120 S 3 S 3 


■ ■ 1 sj 

5040 1 


T 8,l “ ^8 


8,2 


8,3 - 2 


8,4 


8,5 


8,7 


*1 S 7 + S 2 S 6 + S 3 5 S + \ S 4 

Is 2 s 6 + S lS2 S 5 + S lS3 S 4 + I s 2 s 4 + |s 2 s 2 
I S 3S 5 + jS 2 S 2 S 4 + Is^s 2 + ^ Si S 2 S 3 + A s< 
i^S 4 + J 8* s 2 s 3 + £s»8| 

T55 8 ? *> + * 8 ? 8 ? 


750 8 ? 52 


1 


§8 

lg ' 8 40320 1 


The set of operators T. given here permits one to develop the general perturbations up to the 
eighth order. The extension of the given table and the check computations can be performed with 
the general formulas 
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j -k+ 1 

kT. = V" 1 ST. t , 

J , k / O' j ~cr , k“ 1 

CT= 1 


Taking (15) into account, we can write the operator D as 


D = D., (26) 

j=l 


where we put 


D = X*-^ + Y* • ~ + H*.A 
J J Bx J By J Brj 


In order to abbreviate the writing we introduce the symbols 


= L 


A. cf> - 


k=l 

j-1 

L 

k = l 


T ; _„ k = 1, 2 j-1 


D ,- k ^ j = 2. 3. ... 


0, Aj<^ = 0, 


representing the result of applying operators L. = [T , T 2 , . . . , T 1 , 0, 0, 0 . . . ] , 

A. = [D j, D . . . . Dj, 0, 0, 0 . ■ .] to the decomposition of <t> into series with respect to the 
small parameter. 


From (4) - (6), (15), (22) - (24), (25) and (26), we deduce the set of the partial differential 


equations 

Sa . 

A -i = X. - X* + L.X - A. a , 

By j j J J 

(27) 


Bb. 

A- — = Y. - Y* + T.A + LY - A.b, 

dy ■> J J J J 

(28) 
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B|5. 

A = H. - H* + L. H - A. p. 

By J J J J 


For the effects of the first order we have 


MX, 


Tl A 


MY, 


MH, 


A 


Bo i 

By 




PH, 


( 29 ) 


and 


a, = QPX r bj = QPYj, Pj = OPHj . 

Taking into account that the D. operators contain only the long period terms and that a . , b. , p. 
contain only the short period terms, we conclude that A. a , A. b , A. p contain no long period terms 


and thus, to avoid the secular terms in a , b , p , we have to put 

X* = ML. x, (30) 

J J 

Y* = T.A + ML. Y, (31) 

H* = ML. H (32) 

It follows from (27) - (29) and (30) - (32): 

□ . = 0P(X. + L x _ A a), (33) 

J J J J 

b j * QP (Y. + L. Y - A. b), (34) 

p. = 0P(H. + L. H - A. P). (35) 


Evidently (30) - (35) answer the question as to how the long period terms will be formed in higher 
approximations in the Krylov- Bogoliubov method either directly or as a result of the "cross- 
action" of the short period effects. 

These equations can be written in a somewhat simpler form if the series for X, Y and H are 
reduced to one term only. Then we have L. $ = T 0 , and the basic equations become 

X* = MT. , X, 
j i i 
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mu 


ilium mi 


Y* = T.A + MT. , Y, 

j i j i 

H* = MT. , H, 
j j i 


and 


a . = QP(T. , X - A.a), 

j j 1 J 

b. = QP^Y - A.b), 
p. = QPCT.^H - A. p). 

HIGHER ORDER PERTURBATIONS IN THE POINCAR^ AND VON ZEIPEL METHOD 

The introduction of the partial differential operators T and L. permits us also to write the 
equations of the Poincare and von Zeipel method in a very concise form. Consider the system of 
canonical equations: 


dx 

3F 

dy _ 

BF 

d t 

+ V 

dt 

Bx 

d£ 

BF 

_ 

BF 

d t 

3t) ’ 

d t 



We assume that the Hamiltonian F is developable in powers of a small parameter and has the form 


F = F 0 (x) + Fj(x, §; y, n) + F 2 (x, § ; y, n) + . . . 


(36) 


The functions F. (j = 1, 2, 3, • • • ) are periodic in y and n with period 2n in each component: 

F j,n,vfr’5') expi(n-y + v-t])- 


F. = 

J 


We shall determine a canonical transformation 


3S 

. _ 3S 

’ 

y — y + 1 

Bx* 

3S 

3S 

’ 

1* = 1 + — 
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S(x*, g*; y, S j ^ x *’ 5 * : y ’ 

i=i 


such that the new Hamiltonian 


F * = F 0 + F 1 + F 2 + 


(37) 


does not contain the short period argument y* ; in other words, such that the condition 

F(x, g; y, n ) = F * (**• §* : -■ '!*'> 

is satisfied. 

Putting 


(38) 


h(x*, g*; y, r,) = 


dS 


s as 

k (x*, g*; y, Tp = — — ■ 
ox 


x(x*. g*; y, n) = 

on 


3S 

1 


«(x*. g*; y, n) = 


as 

dp 


we write (38) as 

F [x* + h(x*, §*; y, ni.5* + x<v. 5*; y- nh y. ’ll = F * [x*, g*; n + k(x*. g*; y, n i] 


or, changing the notation, as 

f [x + h(x, g, y, tii. g + x (*’ 5; y. i); y. n] = F * [*> -- n + •<(*- 5; y, n)] • 


(39) 


Introducing the displacement operators 

1 + T(x , g; y, n) = exp 

and 


•> (x, g; y. nl • + x(*> S; 

Ox 


■ 5i " "''#]• 


1 + T*(x, g; y, ti) = exp k (x , g; y, t|) 
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we can write (39) in the form 


[l + T( x , §; y, n)] F(x, §; y, t)) = [l + T* (x , §; y, ti)] F*(x, f; n). 


(40) 


Let us define the operators s k and s* by means of the equations 


3S k B 

k By Bx + Bt| B§ 

BS. 3 

o* __ k o 

k B^ Bi| 

Then we have, similarly as before, 


CO 



3 = 1 



(41) 


€0 CD 


1 + T* = exp y ' 
3 = 1 

S* = ^ T*, 

3 = 0 

(42) 

L j = [V.-Vi- 

, Tj. 0, 0, . . .] 


L * = (5-3.T 

. TJ, 0, 0, • • •] • 



The operators T. are expressible in terms of S k , and the operators T* in terms of 6* , by means 
of the formulas given in the previous section. Making use of (37), (38), (41), and (42), we obtain 


L 

j =0 k=0 


£ w. ■ L 


j=0 k=0 


3 

y 1 t* f*. 

/ , k ’ 


K = F 0 (x), 


T* F* = 0, j > 0, 


and consequently 


V -1 


( J 

k=0 


Vk F k = 


y-. r j = i, 2, 3, • • • 

E r >-- K i 


(43) 
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Taking into account 


we can rewrite (43) as 


T j F 0 = 


By Bx 


cr. 


s i^ Fo- 


bs. 

A . — - + 3>. = F* , 

By J J 


(44) 


where 


A 


BF„ 


Bx 


and 


<J>. = F. + (T. - 8.) F„ + L.F - L*F*. 
j i ' i i o T j j 

From the system of the linear partial differential equations (44) we can determine the func- 
tions S. and F* step by step. In order to dispose of the secular terms in S. we have to put F* = 
M ; then we obtain 

j ’ 


S. = OP-1).. 

J J 

The system of transformed equations becomes 

(45) 

(46) 


dx* 

BF* 

= o, 

dy* _ 

d t 

= + v* 7 


d t 


BF* 

Bn*"' 

d n* 

i 

i 

Qj 

* 

d t 

d t 



Besides the integral of energy, the new system also possesses the integral x* = const. 

The system (46) can be integrated independently from the system (45) and after the integration 
the angle y* can be obtained by a plain quadrature. 

A further reduction is possible if F* can be re-arranged so that the purely secular term is of 
lower order than the periodic terms. By repeating the process of elimination of the periodic 
terms, we can obtain the solution of the original problem in the form of a Fourier series with 
arguments linear with respect to time. Such a reduction is possible in the case of the artificial 
satellite of the earth (Brouwer, 1959), but it is not always possible in the lunar problem or in the 
stellar three body problem. If the close companion (the lunar orbiter) is in a highly eccentric orbit 
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and the osculating plane has a high inclination toward the orbital plane of the distant companion 
(Brown, 1936, 1937), then the solution in the form of standard trigonometric series generally 
speaking cannot be obtained. 


CONCLUSION 

The method of Krylov and Bogoliubov 'does not presuppose that the forces are conservative. 
Thus, the importance and the generality of this method are quite evident. The system of the dif- 
ferential operators and the algorithm given here permit the computation of the higher order effects 
up to any order. The process is formal and from the standpoint of pure mathematics might suffer, 
as all astronomical theories do, from the presence of small divisors. 

Recently, the method of Krylov and Bogoliubov was successfully applied by Struble (1961) and 
by Kyner (1965) to the problem of motion of the artificial satellite. The author of the present paper 
has applied it to the problem of the motion of a lunar orbiter, and an exposition of the results will 
appear in a later paper. 
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